Standard model for liquid water withstands x-ray probe 
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We present a series of ab-initio calculations of spectroscopic properties of liquid water at ambient 
conditions. Our results show that all available theoretical and experimental evidence is consistent 
with the standard model of the liquid as comprising molecules with approximately four hydrogen 
bonds. In particular, this model cannot be discounted on the basis of comparisons between measured 
and computed x-ray absorption spectra (XAS), as recently suggested. Our simulations of ice XAS 
including the lowest lying excitonic state are in excellent agreement with experiment and those of 
the TIP4P model of water are in reasonable agreement with recent measurements. Hence we propose 
that the standard, quasi-tetrahedral model of water, although approximate, represents a reasonably 
accurate description of the local structure of the liquid. 

PACS numbers: 61.10.Ht,78.70.Dm,61.20. Ja,61.20.Gy 
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Understanding the hydrogen-bonding in liquid water 
is fundamental to a thorough comprehension of the driv- 
ing forces behind many physical, chemical and biologi- 
cal processes. For the past forty years, water has been 
modeled as having a local structure not dissimilar to 
that of the solid phase. In crystalline ice, each molecule 
is hydrogen-bonded to exactly four others, in an ap- 
proximately tetrahedral arrangement i£ In the standard 
model of the liquid the local, molecular coordination is 
also approximately tetrahedral, with about 3.6 hydrogen 
bonds/molecule^ Classical potentials, fitted to x-ray and 
neutron scattering data and to measured thermodynamic 
properties, yield a quasi-tetrahedral local structure of the 
liquid when used in molecular dynamics (MD) simula- 
tions, and simultaneously reproduce many other phys- 
ical properties X More recently, first principles, density 
functional theory*^ (DFT) calculations have confirmed 
the quasi tetrahedral coordination of water molecules in 
the liquid. However, when neglecting proton quantum 
effects, converged DFT-MD simulations produce a more 
ice-like structure for the liquid^ than that inferred from 
experimental estimates of the radial distribution func- 
tions. 3 

X-ray absorption spectroscopy techniques may help to 
elucidate the local structure of liquid water. In these 
approaches, high-energy, x-ray photons excite electrons 
from deep atomic core levels to states above the Fermi 
level. If structural models are available, together with 
their electronic structure, it is possible to simulate the x- 
ray absorption spectra (XAS) of several candidate struc- 
tures and establish which one best represent the mea- 
surements. The success of this approach hinges on the 
accuracy of XAS simulations and that is the focus of this 
letter. 

Recent XAS experiments on water have questioned the 
standard, tetrahedral model of the liquid. 9 These experi- 
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ments have compared the XAS of bulk ice, the ice surface 
and water and concluded that the liquid contains sig- 
nificantly more broken hydrogen-bonds than previously 
thought. The spectra are qualitatively consistent with 
independent experiments^^ however, their interpreta- 
tion remains controversial. In conjunction with DFT cal- 
culations of the XAS, the experimental results have led 
to the proposal of a new structural model in which each 
molecule is strongly hydrogen-bonded to only two oth- 
ers. These conclusions have major implications for the 
physics and chemistry of liquid water. 

In this work, we report on first principles DFT sim- 
ulations of the XAS of ice and liquid water at ambient 
pressures. Our approach can reliably predict the near K 
edge stucture of oxygen in crystalline ice. Given the very 
good agreement between theory and experiment found in 
the case of ice, we used the same approach to simulate 
the XAS of the liquid as described by the TIP4P classical 
potential;- this potential yields a local structure consis- 
tent with the standard model. Our results are in reason- 
able agreement with experimental XAS. Further analysis 
of molecular species of the liquid with broken hydrogen- 
bonds shows that their spectral signature is qualitatively 
different from the average XAS of the liquid; hence in- 
creasing the proportion of these species would yield spec- 
tra which are significantly different from the measured 
ones. We also discuss previous theoretical approaches to 
simulate the XAS of water and we show that our re- 
sults and those present in the literature are all consistent 
with a quasi-tetrahedral model of the local structure of 
the liquid. 

The x-ray absorption cross section is calculated to first 
order using Fermi's golden rule: 

a(u) = 4ir 2 a hujY^ \M t ^ f \ 2 6(E f - E, - %uS) , (1) 
/ 

where huj is the energy of the absorbed photon, which 
should match the energy difference Ef — Ei between ini- 
tial and final electronic states; cto is the fine structure 
constant; and Mi_>/ are the matrix elements of the tran- 
sition between initial and final states: and \^f) : 
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evaluated within the electric-dipole approximation as 

M^ f = |e ■ R| %) » 5 % |e • r| ^} , (2) 

where e is the polarization direction of the electromag- 
netic vector potential, R and r are the many-electron and 
single-electron position operators respectively, and 
refer to the pair of single-particle states involved in the 
transition. Here the initial state is fixed as the Is eigen- 
state of the oxygen atom. Single-particle approximations 
of the many-electron matrix elements are accurate up to 
a factor S, approximately constant for all transitions. 13 

In our approach the final state of the electronic sys- 
tem is calculated in the presence of the core hole which 
results from the x-ray excitation. To reduce the compu- 
tational cost, we occupy only the first available empty 
band with the excited electron and relax the DFT elec- 
tronic structure within this constraint. 14 We then use the 
corresponding self-consistent potential to generate unoc- 
cupied levels higher in energy than the first excitation. 
We use the pseudopotential approximation, and model 
the x-ray excited atom with a pseudopotential derived 
from an oxygen atom with one electron removed from the 
Is level. We employ norm-conserving Hamann pseudopo- 
tentials 1 ^ for oxygen and a Gaussian pseudopotential for 
hydrogen by Giannozzi. The breaking of spin-degeneracy 
accompanying a single-electron excitation is not included. 
The matrix elements between the atomic core 

level and the excited conduction band are computed us- 
ing a frozen-core approximationii2ii£ii£ We employ DFT 
within the generalized-gradient-approximation using the 
Perdew-Burke-Ernzerhof exchange-correlation functional 
(PBEGGA) 1 ^ and the PWSCF code. 19 . The single particle 
wavefunctions are expanded in a plane-wave basis with 
an energy cut-off of 85 Ryd. 

Recent NEXAFS experiments^*^ reported x-ray ab- 
sorption of crystalline ice samples prepared by molecular 
deposition on the Pt[lll] surface. Helium scattering in- 
dicates^ 1 - the presence of molecular crystals with surfaces 
consistent with either the (0001) surface of hexagonal ice 
(Ih) or the (111) surface of cubic ice (Ic). (Recent work 22 
indicates that coverage of crystalline ice on Pt[lll] may 
not be uniform) . In Fig^we compare the NEXAFS mea- 
surements with our calculated XAS. 

We used lattice constants of ice Ih yielding a density 
of 1.00 g cm~ 3 with a c/a ratio of 0.945. The density 
of ice Ic was fixed at the experimental value of 0.931 g 
cnr 3 using c/a = l.0m& The primitive cell atomic struc- 
tures were relaxed, with these constraints on the unit cell 
volumes, until the forces computed using PBEGGA were 
less than 10~ 3 a.u. 

The primitive cell of ice Ih (Ic) contains 12 (8) wa- 
ter molecules and we constructed a 2 x 2 x 2 supercell 
of 96 (64) molecules for our calculations^ The excited 
state electronic structure, including the excited electron 
in the lowest available conduction band, was converged 
using a uniform grid of 27 (64) fc-points in the Brillouin 
zone containing 10 (36) symmetry-inequivalent points. 
All absorption cross sections were averaged over the three 
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FIG. 1: Left: NEXAFS spectrum of crystalline ice I 
from Refp- (blue) and calculated XAS of heaxgonal ice Ih 
(black, solid) and cubic ice Ic (black, dotted). Right: NEX- 
AFS (blue,solid) and XRS (blue,dashed) spectra of water from 
Ref£, and calculated XAS of water (red) from the TIP4P MD 
simulation at 300 K. Vertical thickness of the calculated curve 
is the associated standard error from a sample of 32 water 
molecules. 



Cartesian directions. 

The agreement between calculated and measured XAS 
(see Fig. 1) is excellent, with all qualitative features re- 
produced accurately. The computed spectrum is aligned 
with an onset energy of 535 eV associated with the lowest 
computed transition (at the zone center), broadened us- 
ing a Gaussian lineshape with 0.4 eV standard deviation, 
and renormalized to reproduce the peak height of experi- 
ment. In general the width of the computed peaks is too 
narrow in comparison with experiment. This is consistent 
with the general underestimation of band-width within 
current local approximations to DFT and originates from 
an approximate description of exchange interactions. Im- 
provements are possible using self-interaction corrected 
DFT calculations^ 3 " or by using GW quasi-particle cor- 
rections. 24 Note that ice Ih and Ic differ structurally only 
beyond first nearest neighbor molecules. Consequently, 
their spectra are very similar, given the short spatial 
range of the oxygen K edge transition probability, which 
is proportional to the overlap of the p-character of con- 
duction bands with the oxygen ls-orbital, possessing a 
0.2 A effective radius. 

We model liquid water using the classical, TIP4P point 
charge model^ which qualitatively reproduces the phase 
diagram of water^ and which we have previously used 
to analyze the electronic structure of the liquid^ In this 
model water molecules are treated as rigid. We took 
10 snapshots of 32 water molecules from a 200 ps MD 
simulationSL 2 ^ spaced 20 ps apart. To approximate the 
XAS experiment, we average computed spectra from each 
molecule in a given snapshot and further average over all 
snapshots. In total we combined the results of 320 spec- 
tral calculations to approximate the XAS of the liquid. 
We found very little variation between the XAS of indi- 
vidual, uncorrelated snapshots computed using only the 
zone-center fc-point. Therefore, we picked one represen- 
tative snapshot of 32 molecules, and carried out calcu- 
lations using 27 /c-points in the Brillouin Zone, for each 
excited molecule. Each spectrum of discrete transitions 
is broadened with a Gaussian lineshape with a standard 
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FIG. 2: Left : Comparison of the calculated XAS of water 
(red) and crystalline ice Ih (black). Right: A similar com- 
parison between ice lh and a model of amorphous ice (red). 
Also indicated in both panels are the characteristic spectra 
from water molecules with at least one broken donor hydro- 
gen bond (see text). Vertical thickness indicates standard 
error as in Fig. Q We assume that the experimental data for 
solid and liquid 9 arc quantitatively comparable and use the 
same normalization factor and onset energy for both theoret- 
ical curves. 



deviation of 0.3 eV. The averaged XAS is compared in 
Fig. n with the experimental NEXAFS and XRS spectra 
of water at 300 Km- The vertical width of the theoreti- 
cal curve indicates an estimate of the standard error at 
a given energy. Our calculation is in good qualitative 
agreement with experiment, indicating that the tetrahe- 
dral local model of the liquid accounts reasonably well 
for measured XAS. 

The contribution to the XAS of the liquid from 
molecules with broken or distorted hydrogen-bonds (HB) 
is compared to that of the entire liquid in Fig. Our 
model contains approximately 20% of molecules with bro- 
ken HB, which are defined according to two definitions: 
(I) exceeding a maximum separation of oxygen atoms of 
3.5 A and a maximum angle of the donated hydrogen 
from the line joining both oxygens of 40°; (II) exceeding 
a maximum oxygen-oxygen separation which depends on 
the angle as outlined in Ref£. 

There is a clear, qualitative difference between the spec- 
tra of species with broken HB and the averaged spec- 
trum. Such differences were predicted in Ref£, where it 
was concluded that the liquid should contain about 80% 
of broken HB species in order to produce the measured 
XAS. We do not find this to be necessary when using our 
approach for the calculation of XAS, which differs from 
that used in Ref£ . 

In Fig. |2 we also compare the XAS of ice Ih and of 
a rather primitive model of amorphous or disordered 
ice, generated using the TIP4P potential by quenching 
a 300 K liquid sample of 32 molecules down to 100 K. 
One snapshot at this temperature was used to generate 
the XAS in Fig. [3 The spectra of the broken HB species 
are qualitatively similar in the disordered sample and in 
the liquid. However, only one (two) of the molecules in 
this snapshot had less than two donor HB according to 
definition I (II), and yet the XAS is qualitatively differ- 
ent from that of the crystalline sample. Clearly disorder 
plays a role in reducing the peak height at 539-540 eV. 




FIG. 3: Comparison of DFT XAS within various approxima- 
tions for the x-ray absorption process, viz. excited state core 
hole (XCH), full core hole (FCH) and half core hole (HCH). 
(a) hexagonal ice; (b) amorphous ice; (c) liquid water; (d) 
liquid water sampling only broken-hydrogen-bonded species. 
Note that the calculated HCH spectrum of a standard, quasi- 
tetrahedral liquid water model appears ice-like in compari- 
son with experiment, however it is still clearly distinguishable 
from ice XAS. 



We also observe an increase in intensity around 538 eV 
and at much higher energies above 545 eV. Both of these 
trends are seen in the liquid as well, and can be associ- 
ated with disorder in the oxygen sublattice. On the other 
hand, the significant increases in intensity at the onset 
(535 eV) and the main-edge (537-538 eV) found in the 
liquid can be associated with broken HB, of which there 
are approximately 4-6 times as many in the liquid sam- 
ple. There are therefore two important effects accounting 
for the difference between XAS of ice and XAS of the liq- 
uid: broken HB and disorder of the oxygen sublattice. 

In the recent literature, two approaches have been 
applied to the simulation of the XAS of water in vari- 
ous phases: the full core hole (FCH) and half-core hole 
(HCH) approaches. The FCH technique — is analogous 
to the one adopted here, except that there is no self- 
consistent inclusion of the excited electron. Such an ap- 
proach should accurately reproduce the higher energy ex- 
citations to progressively more delocalizcd states, and at 
high energy yield the same results as our approach, which 
we call excited-state-core- hole (XCH). However, account- 
ing for the impact of a localized excited electron at low 
energies leads to differences between FCH and XCH re- 
sults; these are shown in Fig. FCH tends to overesti- 
mate the absorption cross section at and near the onset, 
consequently underestimating the main peak height, by 
virtue of the oscillator strength sum-rule. 

In the HCH approach only half of an electron is re- 
moved from the core-state, thereby simulating a transi- 
tion state in the x-ray excitation process. We have found 
that in the presence of the HCH pseudopotential, inclu- 
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sion of an excited electron at the conduction band mini- 
mum has almost no impact on the calculated XAS. The 
HCH spectra tend to overemphasize the main peak and 
underestimate the excitonic near-edge intensity. This is 
consistent with the reduced binding energy of the half- 
core-hole in the excited oxygen atom. 

We note that the same spectral trends exist in all ap- 
proaches (FCH, XCH and HCH) with respect to increas- 
ing disorder and number of broken HB: a decrease in the 
main peak height and an increase in intensity at the near- 
edge are observed. However, in the limit of converged 
fc-point sampling and when employing homogeneous nu- 
merical broadenings, we find that only the XCH approach 
provides a consistent agreement with experiment. 

In summary, using DFT electronic structure calcula- 
tions, we find excellent agreement between the calcu- 
lated and experimentally measured XAS of ice. Using 
the TIP4P classical potential to simulate a standard, 
quasi-tetrahedral model of liquid water we find reason- 
able agreement between our DFT calculations and exper- 
iment. Furthermore, restriction of our analysis to those 



molecular species with broken HB leads to simulated 
XAS which are qualitatively different from experiment. 
Therefore we conclude that the percentage of broken HB 
in the standard model is consistent with experimental 
observations. Finally, despite the different spectra pro- 
duced by various approximations used to simulate the 
x-ray absorption process, we find no conclusive evidence 
within each approach to discount the standard model of 
the liquid. 
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